/* glsl.h - vector and matrix classes to mimic the functionality
 *   of the opengl shading language math objects in C++.
 *
 * Copyright (C) 2004 Bert van der Weerd <bsod@well.com>. Licenced
 * under the Lesser GNU Public Licence (LGPL), version 2 or later. 
 *
 *
 * You really want to use all possible optimization, especially 
 * inlining, and SSE2 code generation on x86 for these classes 
 * to be implemented efficiently. Everyting is defined as inline
 * templates, so you can choose when optimizing to go for code bloat,
 * which is often faster, or to have everything as functions.
 * 
 * Unlike it's name suggests, this code does not depend on or
 * use OpenGL in any way.
 *
 * There are things in the GL Shading Language that cannot be
 * expressed in C++. These are the SWIZZLE OPERATIONS and COMPONENT
 * NAME ALIASSES.
 *
 * The only functions NOT implemented are the texture functions, the
 * fragment shader specific functions, and the noise functions.
 *  
 * Basically this file defines the following classes:
 *
 * vec2 - vector of two floats; v[2] (x,y)
 * vec3 - vector of three floats; v[3] (x,y,z)
 * vec4 - vector of four floats; v[4] (x,y,z,w)
 *
 * ivec2, ivec3, ivec4 - vectors of ints
 * bvec2, bvec3, bvec4 - vectors of booleans
 *
 * mat2 - 2x2 matrix of floats
 * mat3 - 3x3 matrix of floats
 * mat4 - 4x4 matrix of floats
 *
 * All these classes have the usual operations defined on them, and
 * operate usually component-wise. vec3(1)/vec3(2) yields a vec3 with
 * all elements 0.5, for example. The vector multiply operator
 * is component-wise also. For dot product use the dot() function.
 *
 * The major difference between this implementation and the real thing
 * is that it uses templates. vec3, ivec3 and
 * bvec3 are all instantiations of the same template vec3_t<T>, so it
 * is possible to write things like bvec3 t = normalize(ceil(bvec(1)))
 * which are obviously not allowed in the standard. You will get a
 * lot of complaints from the compiler tho, giving you a hint that this
 * is not the way to do it. But it also allows you to use
 * vec3_t<double> or even vec3_t<complex<T> >, which could be useful
 * sometimes. The template approach also reduces the amount of typing
 * for me, a good reason for templates indeed.
 *
 * The not() builtin function isn't defined when using gnu c++ because
 * it complains about 'not'  being a reserved word or something.
 *
 * The matrix classes are implemented row-major order, using a single
 * index, as usual in opengl.
 *
 * Additional functions:
 *   transpose() - transpose a matrix
 *   diag()      - make diagonal matrix
 *   scale()     - make a matrix that implements glScale(). diag() and
 *                 scale() are equivalent
 */

#ifndef INCLUDED_GLSL_H
#define INCLUDED_GLSL_H

#ifndef __cplusplus
# error Who lives under the sea in a pineapple?
#endif

#include <cmath>

// we still need to define M_PI and friends
#ifndef M_PI
# define M_PI 3.14159265358979323846
#endif
#ifndef M_LN2
# define M_LN2 0.693147180559945309417
#endif


// GNU C #defines log2, we dont want that
#ifdef __GNUC__
# ifdef log2
#  undef log2
# endif
#endif

// if min and max are #defines, remove them
#ifdef min
#undef min
#endif
#ifdef max
#undef max
#endif


// forward declaration of types for conversion constructors
template <typename T> struct vec2_t;
template <typename T> struct vec3_t;
template <typename T> struct vec4_t;

/*
 * built-in functions for scalars
 */


inline float radians(float degrees) { return static_cast<float>(M_PI/180)*degrees; }
inline float degrees(float radians) { return static_cast<float>(180/M_PI)*radians; }
inline float exp2(float x) { return static_cast<float>(exp(x * M_LN2)); }
inline float log2(float x) { return static_cast<float>(log(x) / M_LN2); }
inline float inversesqrt(float x) { return 1/sqrt(x); }
inline float sign(float x) { return (x>0) ? 1.0f : ((x<0) ? -1.0f:0.0f); }
inline float fract(float x) { return x - floor(x); }
inline float mod(float x, float y) { return fmod(x, y); }
inline float min(float x, float y) { return y < x ? y : x; }
inline float max(float x, float y) { return x < y ? y : x; }
inline float clamp(float x, float minVal, float maxVal) { return min(max(x,minVal),maxVal); }
inline float mix(float x, float y, float a) { return x*(1-a) + y*a; }
inline float step(float edge, float x) { return x<=edge ? 0.0f : 1.0f; }
inline float smoothstep(float edge0,float edge1,float x) { float t = clamp((x-edge0) / (edge1-edge0), 0, 1); return t * t * (3 - 2*t); }
inline float length(float x) { return sqrtf(x*x); }
inline float distance(float p0, float p1) { return length(p0-p1); }
inline float dot(float x, float y) { return x*y; }
inline float normalize(float) { return 1; } // this is not the most useful function in the world
inline float faceforward(float N, float l, float Nref) { return Nref*l < 0 ? N : -N; }
inline float reflect(float I, float N) { return I - 2*N*I*N; }


/*
 * vec2_t<T> - vec2, ivec2, bvec2
 */


template <typename T>
   struct vec2_t
{
	T v[2];

	inline vec2_t(){ v[0]=v[1]=0; }
	inline vec2_t(T t) { v[0] = v[1] = t; }
	inline vec2_t(T x, T y) {v[0]=x; v[1]=y;}
	inline vec2_t(vec3_t<T> v);

	inline T& operator[](int i) { return v[i]; }
	inline T& operator[](int i) const { return v[i]; }
};
typedef struct vec2_t<float> vec2;
typedef struct vec2_t<int>   ivec2;
typedef struct vec2_t<bool>  bvec2;

// operators
template <typename T>
   inline vec2_t<T> operator+(vec2_t<T> a, vec2_t<T> b) { return vec2_t<T>(a.v[0]+b.v[0], a.v[1]+b.v[1]); }
template <typename T>
   inline vec2_t<T> operator-(vec2_t<T> a, vec2_t<T> b) { return vec2_t<T>(a.v[0]-b.v[0], a.v[1]-b.v[1]); }
template <typename T> // componentwise multiply, use dot() for inner product
   inline vec2_t<T> operator*(vec2_t<T> a, vec2_t<T> b) { return vec2_t<T>(a.v[0]*b.v[0],a.v[1]*b.v[1]); }
template <typename T>
   inline vec2_t<T> operator/(vec2_t<T> a, vec2_t<T> b) { return vec2_t<T>(a.v[0]/b.v[0],a.v[1]/b.v[1]); }

template <typename T>
   inline vec2_t<T> operator*(T f,vec2_t<T> v) { return vec2_t<T>(f*v.v[0], f*v.v[1]); }
template <typename T>
   inline vec2_t<T> operator*(vec2_t<T> v,T f) { return vec2_t<T>(f*v.v[0], f*v.v[1]); }
template <typename T>
   inline vec2_t<T> operator+(T f,vec2_t<T> v) { return vec2_t<T>(f+v.v[0], f+v.v[1]); }
template <typename T>
   inline vec2_t<T> operator+(vec2_t<T> v,T f) { return vec2_t<T>(f+v.v[0], f+v.v[1]); }
template <typename T>
   inline vec2_t<T> operator/(vec2_t<T> v,T f) { return vec2_t<T>(v.v[0]/f, v.v[1]/f); }

template <typename T>
   inline vec2_t<T>& operator+=(vec2_t<T>& a, vec2_t<T> b) { a = a+b; return a; }
template <typename T>
   inline vec2_t<T>& operator-=(vec2_t<T>& a, vec2_t<T> b) { a = a-b; return a; }
template <typename T>
   inline vec2_t<T>& operator*=(vec2_t<T>& a, vec2_t<T> b) { a = a*b; return a; }
template <typename T>
   inline vec2_t<T>& operator/=(vec2_t<T>& a, vec2_t<T> b) { a = a/b; return a; }
template <typename T>
   inline vec2_t<T>& operator+=(vec2_t<T>& a, T f) { a = a+f; return a; }
template <typename T>
   inline vec2_t<T>& operator-=(vec2_t<T>& a, T f) { a = a-f; return a; }
template <typename T>
   inline vec2_t<T>& operator*=(vec2_t<T>& a, T f) { a = a*f; return a; }
template <typename T>
   inline vec2_t<T>& operator/=(vec2_t<T>& a, T f) { a = a/f; return a; }

template <typename T>
   inline vec2_t<T> operator-(vec2_t<T> v) { return vec2_t<T>(-v.v[0],-v.v[1]); }

// built-in functions
template <typename T>
   inline vec2_t<T> radians(vec2_t<T> v) { return vec2_t<T>(radians(v.v[0]),radians(v.v[1])); }
template <typename T>
   inline vec2_t<T> degrees(vec2_t<T> v) { return vec2_t<T>(degrees(v.v[0]),degrees(v.v[1])); }
template <typename T>
   inline vec2_t<T> sin(vec2_t<T> v) { return vec2_t<T>(sin(v.v[0]),sin(v.v[1])); }
template <typename T>
   inline vec2_t<T> cos(vec2_t<T> v) { return vec2_t<T>(cos(v.v[0]),cos(v.v[1])); }
template <typename T>
   inline vec2_t<T> tan(vec2_t<T> v) { return vec2_t<T>(tan(v.v[0]),tan(v.v[1])); }
template <typename T>
   inline vec2_t<T> asin(vec2_t<T> v) { return vec2_t<T>(asin(v.v[0]),asin(v.v[1])); }
template <typename T>
   inline vec2_t<T> acos(vec2_t<T> v) { return vec2_t<T>(acos(v.v[0]),acos(v.v[1])); }
template <typename T>
   inline vec2_t<T> atan(vec2_t<T> v) { return vec2_t<T>(atan(v.v[0]),atan(v.v[1])); }
template <typename T>
   inline vec2_t<T> atan(vec2_t<T> y, vec2_t<T> x) { return vec2_t<T>(atan2(y.v[0],x.v[0]),atan2(y.v[1],x.v[1])); }

template <typename T>
   inline vec2_t<T> pow(vec2_t<T> a, vec2_t<T> b) { return vec2_t<T>(pow(a.v[0],b.v[0]),pow(a.v[1],b.v[1])); }
template <typename T>
   inline vec2_t<T> exp2(vec2_t<T> v) { return vec2_t<T>(exp2(v.v[0]),exp2(v.v[1])); }
template <typename T>
   inline vec2_t<T> log2(vec2_t<T> v) { return vec2_t<T>(log2(v.v[0]),log2(v.v[1])); }
template <typename T>
   inline vec2_t<T> sqrt(vec2_t<T> v) { return vec2_t<T>(sqrt(v.v[0]),sqrt(v.v[1])); }
template <typename T>
   inline vec2_t<T> inversesqrt(vec2_t<T> v) { return vec2_t<T>(inversesqrt(v.v[0]),inversesqrt(v.v[1])); }

template <typename T>
   inline vec2_t<T> abs(vec2_t<T> v) { return vec2_t<T>(abs(v.v[0]),abs(v.v[1])); }
template <typename T>
   inline vec2_t<T> sign(vec2_t<T> v) { return vec2_t<T>(sign(v.v[0]),sign(v.v[1])); }
template <typename T>
   inline vec2_t<T> floor(vec2_t<T> v) { return vec2_t<T>(floor(v.v[0]),floor(v.v[1])); }
template <typename T>
   inline vec2_t<T> ceil(vec2_t<T> v) { return vec2_t<T>(ceil(v.v[0]),ceil(v.v[1])); }
template <typename T>
   inline vec2_t<T> fract(vec2_t<T> v) { return vec2_t<T>(fract(v.v[0]),fract(v.v[1])); }
template <typename T>
   inline vec2_t<T> mod(vec2_t<T> v, T f) { return vec2_t<T>(mod(v.v[0],f),mod(v.v[1],f)); }
template <typename T>
   inline vec2_t<T> mod(vec2_t<T> a, vec2_t<T> b) { return vec2_t<T>(mod(a.v[0],b.v[0]),mod(a.v[1],b.v[1])); }
template <typename T>
   inline vec2_t<T> min(vec2_t<T> v, vec2_t<T> u) { return vec2_t<T>(min(v.v[0],u.v[0]),min(v.v[1],u.v[1])); }
template <typename T>
   inline vec2_t<T> min(vec2_t<T> v, T f) { return vec2_t<T>(min(v.v[0],f),min(v.v[1],f)); }
template <typename T>
   inline vec2_t<T> clamp(vec2_t<T> v, T minVal, T maxVal) { return vec2_t<T>(clamp(v.v[0],minVal,maxVal),clamp(v.v[1],minVal,maxVal)); }
template <typename T>
   inline vec2_t<T> clamp(vec2_t<T> v, vec2_t<T> minVal, vec2_t<T> maxVal) {
	return vec2_t<T>(clamp(v.v[0],minVal.v[0],maxVal.v[0]),clamp(v.v[1],minVal.v[1],maxVal.v[1]));
}
template <typename T>
   inline vec2_t<T> mix(vec2_t<T> x, vec2_t<T> y, float a) { return vec2_t<T>(mix(x.v[0],y.v[0],a),mix(x.v[1],y.v[1],a)); }
template <typename T>
   inline vec2_t<T> mix(vec2_t<T> x, vec2_t<T> y, vec2_t<T> a) { return vec2_t<T>(mix(x.v[0],y.v[0],a.v[0]),mix(x.v[1],y.v[1],a.v[1])); }
template <typename T>
   inline vec2_t<T> step(vec2_t<T> edge, vec2_t<T> x) { return vec2_t<T>(step(edge.v[0],x.v[0]),step(edge.v[1],x.v[1])); }
template <typename T>
   inline vec2_t<T> smoothstep(vec2_t<T> edge0,vec2_t<T> edge1,vec2_t<T> x) {
	return vec2_t<T>(smoothstep(edge0.v[0],edge1.v[0],x.v[0]),smoothstep(edge0.v[1],edge1.v[1],x.v[1]));
}

template <typename T>
   inline T length(vec2_t<T> v) { return static_cast<T>(sqrtf(static_cast<float>(dot(v,v)))); }
template <typename T>
   inline T distance(vec2_t<T> p0, vec2_t<T> p1) { return length(p0-p1); }
template <typename T>
   inline T dot(vec2_t<T> a, vec2_t<T> b) { return a.v[0]*b.v[0] + a.v[1]*b.v[1]; }
template <typename T>
   inline vec2_t<T> normalize(vec2_t<T> v) { T r = length(v); return vec2_t<T>(v.v[0]/r, v.v[1]/r); }
template <typename T>
   inline vec2_t<T> faceforward(vec2_t<T> N, vec2_t<T> I, vec2_t<T> Nref) {	return dot(Nref,I) < 0 ? N : -N; }
template <typename T>
   inline vec2_t<T> reflect(vec2_t<T> I, vec2_t<T> N) { return I - 2*dot(N,I)*N; }

template <typename T>
   inline bvec2 lessThan(vec2_t<T> a, vec2_t<T> b) { return bvec2(a.v[0]<b.v[0],a.v[1]<b.v[1]); }
template <typename T>
   inline bvec2 lessThanEqual(vec2_t<T> a, vec2_t<T> b) { return bvec2(a.v[0]<=b.v[0],a.v[1]<=b.v[1]); }
template <typename T>
   inline bvec2 greaterThan(vec2_t<T> a, vec2_t<T> b) { return bvec2(a.v[0]>b.v[0],a.v[1]>b.v[1]); }
template <typename T>
   inline bvec2 greaterThanEqual(vec2_t<T> a, vec2_t<T> b) { return bvec2(a.v[0]>=b.v[0],a.v[1]>=b.v[1]); }
template <typename T>
   inline bvec2 equal(vec2_t<T> a, vec2_t<T> b) { return bvec2(a.v[0]==b.v[2],a.v[1]==b.v[1]); }
template <typename T>
   inline bvec2 notEqual(vec2_t<T> a, vec2_t<T> b) { return bvec2(a.v[0]!=b.v[2],a.v[1]!=b.v[1]); }

inline bool any(bvec2 v) { return v.v[0] || v.v[1]; }
inline bool all(bvec2 v) { return v.v[0] && v.v[1]; }
#ifndef __GNUC__
inline bvec2 not(bvec2 v) { return bvec2(!v.v[0],!v.v[1]); }
#endif


/*
 * vec3_t<T> - 3vector: vec3, ivec3, bvec3
 */


template <typename T>
   struct vec3_t
{
	T v[3];

	inline vec3_t() { v[0]=v[1]=v[2]=0; }
	inline vec3_t(T t) { v[0]=v[1]=v[2]=t; }
	inline vec3_t(vec2_t<T> v, T z) {v[0]=v.v[0]; v[1]=v.v[1]; v[2] = z; } // vec2 extension
	inline vec3_t(T x, T y, T z) {v[0]=x;v[1]=y;v[2]=z;}
	inline vec3_t(vec4_t<T> v);

	inline T& operator[](int i) { return v[i]; }
	inline T& operator[](int i) const { return v[i]; }
};

typedef struct vec3_t<float> vec3;
typedef struct vec3_t<int>   ivec3;
typedef struct vec3_t<bool>  bvec3;

// operators
template <typename T>
   inline vec3_t<T> operator+(vec3_t<T> a, vec3_t<T> b) { return vec3_t<T>(a.v[0]+b.v[0], a.v[1]+b.v[1], a.v[2]+b.v[2]); }
template <typename T>
   inline vec3_t<T> operator-(vec3_t<T> a, vec3_t<T> b) { return vec3_t<T>(a.v[0]-b.v[0], a.v[1]-b.v[1], a.v[2]-b.v[2]); }
template <typename T> // componentwise multiply, use dot() for inner product
   inline vec3_t<T> operator*(vec3_t<T> a, vec3_t<T> b) { return vec3_t<T>(a.v[0]*b.v[0],a.v[1]*b.v[1],a.v[2]*b.v[2]); }
template <typename T>
   inline vec3_t<T> operator/(vec3_t<T> a, vec3_t<T> b) { return vec3_t<T>(a.v[0]/b.v[0],a.v[1]/b.v[1],a.v[2]/b.v[2]); }

template <typename T>
   inline vec3_t<T> operator*(T f,vec3_t<T> v) { return vec3_t<T>(f*v.v[0], f*v.v[1], f*v.v[2]); }
template <typename T>
   inline vec3_t<T> operator*(vec3_t<T> v,T f) { return vec3_t<T>(f*v.v[0], f*v.v[1], f*v.v[2]); }
template <typename T>
   inline vec3_t<T> operator+(T f,vec3_t<T> v) { return vec3_t<T>(f+v.v[0], f+v.v[1], f+v.v[2]); }
template <typename T>
   inline vec3_t<T> operator+(vec3_t<T> v,T f) { return vec3_t<T>(f+v.v[0], f+v.v[1], f+v.v[2]); }
template <typename T>
   inline vec3_t<T> operator/(vec3_t<T> v,T f) { return vec3_t<T>(v.v[0]/f, v.v[1]/f, v.v[2]/f); }

template <typename T>
   inline vec3_t<T>& operator+=(vec3_t<T>& a, vec3_t<T> b) { a = a+b; return a; }
template <typename T>
   inline vec3_t<T>& operator-=(vec3_t<T>& a, vec3_t<T> b) { a = a-b; return a; }
template <typename T>
   inline vec3_t<T>& operator*=(vec3_t<T>& a, vec3_t<T> b) { a = a*b; return a; }
template <typename T>
   inline vec3_t<T>& operator/=(vec3_t<T>& a, vec3_t<T> b) { a = a/b; return a; }
template <typename T>
   inline vec3_t<T>& operator+=(vec3_t<T>& a, T f) { a = a+f; return a; }
template <typename T>
   inline vec3_t<T>& operator-=(vec3_t<T>& a, T f) { a = a-f; return a; }
template <typename T>
   inline vec3_t<T>& operator*=(vec3_t<T>& a, T f) { a = a*f; return a; }
template <typename T>
   inline vec3_t<T>& operator/=(vec3_t<T>& a, T f) { a = a/f; return a; }

template <typename T>
   inline vec3_t<T> operator-(vec3_t<T> v) { return vec3_t<T>(-v.v[0],-v.v[1],-v.v[2]); }

// built-in functions
template <typename T>
   inline vec3_t<T> radians(vec3_t<T> v) { return vec3_t<T>(radians(v.v[0]),radians(v.v[1]),radians(v.v[2])); }
template <typename T>
   inline vec3_t<T> degrees(vec3_t<T> v) { return vec3_t<T>(degrees(v.v[0]),degrees(v.v[1]),degrees(v.v[2])); }
template <typename T>
   inline vec3_t<T> sin(vec3_t<T> v) { return vec3_t<T>(sin(v.v[0]),sin(v.v[1]),sin(v.v[2])); }
template <typename T>
   inline vec3_t<T> cos(vec3_t<T> v) { return vec3_t<T>(cos(v.v[0]),cos(v.v[1]),cos(v.v[2])); }
template <typename T>
   inline vec3_t<T> tan(vec3_t<T> v) { return vec3_t<T>(tan(v.v[0]),tan(v.v[1]),tan(v.v[2])); }
template <typename T>
   inline vec3_t<T> asin(vec3_t<T> v) { return vec3_t<T>(asin(v.v[0]),asin(v.v[1]),asin(v.v[2])); }
template <typename T>
   inline vec3_t<T> acos(vec3_t<T> v) { return vec3_t<T>(acos(v.v[0]),acos(v.v[1]),acos(v.v[2])); }
template <typename T>
   inline vec3_t<T> atan(vec3_t<T> v) { return vec3_t<T>(atan(v.v[0]),atan(v.v[1]),atan(v.v[2])); }
template <typename T>
   inline vec3_t<T> atan(vec3_t<T> y, vec3_t<T> x) { return vec3_t<T>(atan2(y.v[0],x.v[0]),atan2(y.v[1],x.v[1]),atan2(y.v[2],x.v[2])); }

template <typename T>
   inline vec3_t<T> pow(vec3_t<T> a, vec3_t<T> b) { return vec3_t<T>(pow(a.v[0],b.v[0]),pow(a.v[1],b.v[1]),pow(a.v[2],b.v[2])); }
template <typename T>
   inline vec3_t<T> exp2(vec3_t<T> v) { return vec3_t<T>(exp2(v.v[0]),exp2(v.v[1]),exp2(v.v[2])); }
template <typename T>
   inline vec3_t<T> log2(vec3_t<T> v) { return vec3_t<T>(log2(v.v[0]),log2(v.v[1]),log2(v.v[2])); }
template <typename T>
   inline vec3_t<T> sqrt(vec3_t<T> v) { return vec3_t<T>(sqrt(v.v[0]),sqrt(v.v[1]),sqrt(v.v[2])); }
template <typename T>
   inline vec3_t<T> inversesqrt(vec3_t<T> v) { return vec3_t<T>(inversesqrt(v.v[0]),inversesqrt(v.v[1]),inversesqrt(v.v[2])); }

template <typename T>
   inline vec3_t<T> abs(vec3_t<T> v) { return vec3_t<T>(abs(v.v[0]),abs(v.v[1]),abs(v.v[2])); }
template <typename T>
   inline vec3_t<T> sign(vec3_t<T> v) { return vec3_t<T>(sign(v.v[0]),sign(v.v[1]),sign(v.v[2])); }
template <typename T>
   inline vec3_t<T> floor(vec3_t<T> v) { return vec3_t<T>(floor(v.v[0]),floor(v.v[1]),floor(v.v[2])); }
template <typename T>
   inline vec3_t<T> ceil(vec3_t<T> v) { return vec3_t<T>(ceil(v.v[0]),ceil(v.v[1]),ceil(v.v[2])); }
template <typename T>
   inline vec3_t<T> fract(vec3_t<T> v) { return vec3_t<T>(fract(v.v[0]),fract(v.v[1]),fract(v.v[2])); }
template <typename T>
   inline vec3_t<T> mod(vec3_t<T> v, T f) { return vec3_t<T>(mod(v.v[0],f),mod(v.v[1],f),mod(v.v[2],f)); }
template <typename T>
   inline vec3_t<T> mod(vec3_t<T> a, vec3_t<T> b) { return vec3_t<T>(mod(a.v[0],b.v[0]),mod(a.v[1],b.v[1]),mod(a.v[2],b.v[2])); }
template <typename T>
   inline vec3_t<T> min(vec3_t<T> v, vec3_t<T> u) { return vec3_t<T>(min(v.v[0],u.v[0]),min(v.v[1],u.v[1]),min(v.v[2],u.v[2])); }
template <typename T>
   inline vec3_t<T> min(vec3_t<T> v, T f) { return vec3_t<T>(min(v.v[0],f),min(v.v[1],f),min(v.v[2],f)); }
template <typename T>
   inline vec3_t<T> clamp(vec3_t<T> v, T minVal, T maxVal) { return vec3_t<T>(clamp(v.v[0],minVal,maxVal),clamp(v.v[1],minVal,maxVal),clamp(v.v[2],minVal,maxVal)); }
template <typename T>
   inline vec3_t<T> clamp(vec3_t<T> v, vec3_t<T> minVal, vec3_t<T> maxVal) {
	return vec3_t<T>(clamp(v.v[0],minVal.v[0],maxVal.v[0]),clamp(v.v[1],minVal.v[1],maxVal.v[1]),clamp(v.v[2],minVal.v[2],maxVal.v[2]));
}
template <typename T>
   inline vec3_t<T> mix(vec3_t<T> x, vec3_t<T> y, float a) { return vec3_t<T>(mix(x.v[0],y.v[0],a),mix(x.v[1],y.v[1],a),mix(x.v[2],y.v[2],a)); }
template <typename T>
   inline vec3_t<T> mix(vec3_t<T> x, vec3_t<T> y, vec3_t<T> a) { return vec3_t<T>(mix(x.v[0],y.v[0],a.v[0]),mix(x.v[1],y.v[1],a.v[1]),mix(x.v[2],y.v[2],a.v[2])); }
template <typename T>
   inline vec3_t<T> step(vec3_t<T> edge, vec3_t<T> x) { return vec3_t<T>(step(edge.v[0],x.v[0]),step(edge.v[1],x.v[1]),step(edge.v[2],x.v[2])); }
template <typename T>
   inline vec3_t<T> smoothstep(vec3_t<T> edge0,vec3_t<T> edge1,vec3_t<T> x) {
	return vec3_t<T>(smoothstep(edge0.v[0],edge1.v[0],x.v[0]),smoothstep(edge0.v[1],edge1.v[1],x.v[1]),smoothstep(edge0.v[2],edge1.v[2],x.v[2]));
}

template <typename T>
   inline T length(vec3_t<T> v) { return static_cast<T>(sqrtf(static_cast<float>(dot(v,v)))); }
template <typename T>
   inline T distance(vec3_t<T> p0, vec3_t<T> p1) { return length(p0-p1); }
template <typename T>
   inline T dot(vec3_t<T> a, vec3_t<T> b) { return a.v[0]*b.v[0] + a.v[1]*b.v[1] + a.v[2]*b.v[2]; }
template <typename T>
   inline vec3_t<T> cross(vec3_t<T> a,vec3_t<T> b) {
	vec3_t<T> result;
	result.v[0] = a.v[1]*b.v[2] - a.v[2]*b.v[1];
	result.v[1] = a.v[2]*b.v[0] - a.v[0]*b.v[2];
	result.v[2] = a.v[0]*b.v[1] - a.v[1]*b.v[0];
	return result;
}
template <typename T>
   inline vec3_t<T> normalize(vec3_t<T> v) { T r = length(v); return vec3_t<T>(v.v[0]/r, v.v[1]/r, v.v[2]/r); }
template <typename T>
   inline vec3_t<T> faceforward(vec3_t<T> N, vec3_t<T> I, vec3_t<T> Nref) {	return dot(Nref,I) < 0 ? N : -N; }
template <typename T>
   inline vec3_t<T> reflect(vec3_t<T> I, vec3_t<T> N) { return I - 2*dot(N,I)*N; }

template <typename T>
   inline bvec3 lessThan(vec3_t<T> a, vec3_t<T> b) { return bvec3(a.v[0]<b.v[0],a.v[1]<b.v[1],a.v[2]<b.v[2]); }
template <typename T>
   inline bvec3 lessThanEqual(vec3_t<T> a, vec3_t<T> b) { return bvec3(a.v[0]<=b.v[0],a.v[1]<=b.v[1],a.v[2]<=b.v[2]); }
template <typename T>
   inline bvec3 greaterThan(vec3_t<T> a, vec3_t<T> b) { return bvec3(a.v[0]>b.v[0],a.v[1]>b.v[1],a.v[2]>b.v[2]); }
template <typename T>
   inline bvec3 greaterThanEqual(vec3_t<T> a, vec3_t<T> b) { return bvec3(a.v[0]>=b.v[0],a.v[1]>=b.v[1],a.v[2]>=b.v[2]); }
template <typename T>
   inline bvec3 equal(vec3_t<T> a, vec3_t<T> b) { return bvec3(a.v[0]==b.v[2],a.v[1]==b.v[1],a.v[2]==b.v[2]); }
template <typename T>
   inline bvec3 notEqual(vec3_t<T> a, vec3_t<T> b) { return bvec3(a.v[0]!=b.v[2],a.v[1]!=b.v[1],a.v[2]!=b.v[2]); }

inline bool any(bvec3 v) { return v.v[0] || v.v[1] || v.v[2]; }
inline bool all(bvec3 v) { return v.v[0] && v.v[1] && v.v[2]; }
#ifndef __GNUC__
inline bvec3 not(bvec3 v) { return bvec3(!v.v[0],!v.v[1],!v.v[2]); }
#endif


/*
 * vec4_t<T> - 4vector: vec4, ivec4, bvec4
 */


template <typename T>
   struct vec4_t
{
	T v[4];

	inline vec4_t() {v[0]=v[1]=v[2]=v[3]=0;}
	inline vec4_t(T t) {v[0]=v[1]=v[2]=v[3]=t;}
	inline vec4_t(vec3_t<T> v, T w) {v[0]=v.v[0];v[1]=v.v[2];v[2]=v.v[2];v[3]=w;} // vec3 extension
	inline vec4_t(T x, T y, T z, T w) {v[0]=x;v[1]=y;v[2]=z;v[3]=w;}

	inline T& operator[](int i) { return v[i]; }
	inline T& operator[](int i) const { return v[i]; }
};

typedef struct vec4_t<float> vec4;
typedef struct vec4_t<int>   ivec4;
typedef struct vec4_t<bool>  bvec4;

// operators
template <typename T>
   inline vec4_t<T> operator+(vec4_t<T> a, vec4_t<T> b) { return vec4_t<T>(a.v[0]+b.v[0], a.v[1]+b.v[1], a.v[2]+b.v[2], a.v[3]+b.v[3]); }
template <typename T>
   inline vec4_t<T> operator-(vec4_t<T> a, vec4_t<T> b) { return vec4_t<T>(a.v[0]-b.v[0], a.v[1]-b.v[1], a.v[2]-b.v[2], a.v[3]-b.v[3]); }
template <typename T> // componentwise multiply, use dot() for inner product
inline vec4_t<T> operator*(vec4_t<T> a, vec4_t<T> b) { return vec4_t<T>(a.v[0]*b.v[0],a.v[1]*b.v[1],a.v[2]*b.v[2], a.v[3]*b.v[3]); }
template <typename T>
   inline vec4_t<T> operator/(vec4_t<T> a, vec4_t<T> b) { return vec4_t<T>(a.v[0]/b.v[0],a.v[1]/b.v[1],a.v[2]/b.v[2], a.v[3]/b.v[3]); }

template <typename T>
   inline vec4_t<T> operator*(T f,vec4_t<T> v) { return vec4_t<T>(f*v.v[0], f*v.v[1], f*v.v[2], f*v.v[3]); }
template <typename T>
   inline vec4_t<T> operator*(vec4_t<T> v,T f) { return vec4_t<T>(f*v.v[0], f*v.v[1], f*v.v[2], f*v.v[3]); }
template <typename T>
   inline vec4_t<T> operator+(T f,vec4_t<T> v) { return vec4_t<T>(f+v.v[0], f+v.v[1], f+v.v[2], f+v.v[3]); }
template <typename T>
   inline vec4_t<T> operator+(vec4_t<T> v,T f) { return vec4_t<T>(f+v.v[0], f+v.v[1], f+v.v[2], f+v.v[3]); }
template <typename T>
   inline vec4_t<T> operator/(vec4_t<T> v,T f) { return vec4_t<T>(v.v[0]/f, v.v[1]/f, v.v[2]/f, v.v[3]/f); }

template <typename T>
   inline vec4_t<T>& operator+=(vec4_t<T>& a, vec4_t<T> b) { a = a+b; return a; }
template <typename T>
   inline vec4_t<T>& operator-=(vec4_t<T>& a, vec4_t<T> b) { a = a-b; return a; }
template <typename T>
   inline vec4_t<T>& operator*=(vec4_t<T>& a, vec4_t<T> b) { a = a*b; return a; }
template <typename T>
   inline vec4_t<T>& operator/=(vec4_t<T>& a, vec4_t<T> b) { a = a/b; return a; }
template <typename T>
   inline vec4_t<T>& operator+=(vec4_t<T>& a, T f) { a = a+f; return a; }
template <typename T>
   inline vec4_t<T>& operator-=(vec4_t<T>& a, T f) { a = a-f; return a; }
template <typename T>
   inline vec4_t<T>& operator*=(vec4_t<T>& a, T f) { a = a*f; return a; }
template <typename T>
   inline vec4_t<T>& operator/=(vec4_t<T>& a, T f) { a = a/f; return a; }

template <typename T>
   inline vec4_t<T> operator-(vec4_t<T> v) { return vec4_t<T>(-v.v[0],-v.v[1],-v.v[2], -v.v[3]); }

// built-in functions
template <typename T>
   inline vec4_t<T> radians(vec4_t<T> v) { return vec4_t<T>(radians(v.v[0]),radians(v.v[1]),radians(v.v[2]),radians(v.v[3])); }
template <typename T>
   inline vec4_t<T> degrees(vec4_t<T> v) { return vec4_t<T>(degrees(v.v[0]),degrees(v.v[1]),degrees(v.v[2]),degrees(v.v[3])); }
template <typename T>
   inline vec4_t<T> sin(vec4_t<T> v) { return vec4_t<T>(sin(v.v[0]),sin(v.v[1]),sin(v.v[2]),sin(v.v[3])); }
template <typename T>
   inline vec4_t<T> cos(vec4_t<T> v) { return vec4_t<T>(cos(v.v[0]),cos(v.v[1]),cos(v.v[2]),cos(v.v[3])); }
template <typename T>
   inline vec4_t<T> tan(vec4_t<T> v) { return vec4_t<T>(tan(v.v[0]),tan(v.v[1]),tan(v.v[2]),tan(v.v[3])); }
template <typename T>
   inline vec4_t<T> asin(vec4_t<T> v) { return vec4_t<T>(asin(v.v[0]),asin(v.v[1]),asin(v.v[2]),asin(v.v[3])); }
template <typename T>
   inline vec4_t<T> acos(vec4_t<T> v) { return vec4_t<T>(acos(v.v[0]),acos(v.v[1]),acos(v.v[2]),acos(v.v[3])); }
template <typename T>
   inline vec4_t<T> atan(vec4_t<T> v) { return vec4_t<T>(atan(v.v[0]),atan(v.v[1]),atan(v.v[2]),atan(v.v[3])); }
template <typename T>
   inline vec4_t<T> atan(vec4_t<T> y, vec4_t<T> x) { return vec4_t<T>(atan2(y.v[0],x.v[0]),atan2(y.v[1],x.v[1]),atan2(y.v[2],x.v[2]),atan2(y.v[3],x.v[3])); }

template <typename T>
   inline vec4_t<T> pow(vec4_t<T> a, vec4_t<T> b) { return vec4_t<T>(pow(a.v[0],b.v[0]),pow(a.v[1],b.v[1]),pow(a.v[2],b.v[2]),pow(a.v[3],b.v[3])); }
template <typename T>
   inline vec4_t<T> exp2(vec4_t<T> v) { return vec4_t<T>(exp2(v.v[0]),exp2(v.v[1]),exp2(v.v[2]),exp2(v.v[3])); }
template <typename T>
   inline vec4_t<T> log2(vec4_t<T> v) { return vec4_t<T>(log2(v.v[0]),log2(v.v[1]),log2(v.v[2]),lo2(v.v[3])); }
template <typename T>
   inline vec4_t<T> sqrt(vec4_t<T> v) { return vec4_t<T>(sqrt(v.v[0]),sqrt(v.v[1]),sqrt(v.v[2]),sqrt(v.v[3])); }
template <typename T>
   inline vec4_t<T> inversesqrt(vec4_t<T> v) { return vec4_t<T>(inversesqrt(v.v[0]),inversesqrt(v.v[1]),inversesqrt(v.v[2]),inversesqrt(v.v[3])); }

template <typename T>
   inline vec4_t<T> abs(vec4_t<T> v) { return vec4_t<T>(abs(v.v[0]),abs(v.v[1]),abs(v.v[2]),abs(v.v[3])); }
template <typename T>
   inline vec4_t<T> sign(vec4_t<T> v) { return vec4_t<T>(sign(v.v[0]),sign(v.v[1]),sign(v.v[2]),sign(v.v[3])); }
template <typename T>
   inline vec4_t<T> floor(vec4_t<T> v) { return vec4_t<T>(floor(v.v[0]),floor(v.v[1]),floor(v.v[2]),floor(v.v[3])); }
template <typename T>
   inline vec4_t<T> ceil(vec4_t<T> v) { return vec4_t<T>(ceil(v.v[0]),ceil(v.v[1]),ceil(v.v[2]),ceil(v.v[3])); }
template <typename T>
   inline vec4_t<T> fract(vec4_t<T> v) { return vec4_t<T>(fract(v.v[0]),fract(v.v[1]),fract(v.v[2]),fract(v.v[3])); }
template <typename T>
   inline vec4_t<T> mod(vec4_t<T> v, T f) { return vec4_t<T>(mod(v.v[0],f),mod(v.v[1],f),mod(v.v[2],f),mod(v.v[3],f)); }
template <typename T>
   inline vec4_t<T> mod(vec4_t<T> a, vec4_t<T> b) { return vec4_t<T>(mod(a.v[0],b.v[0]),mod(a.v[1],b.v[1]),mod(a.v[2],b.v[2]),mod(a.v[3],b.v[3])); }
template <typename T>
   inline vec4_t<T> min(vec4_t<T> v, vec4_t<T> u) { return vec4_t<T>(min(v.v[0],u.v[0]),min(v.v[1],u.v[1]),min(v.v[2],u.v[2]),min(v.v[3],u.v[3])); }
template <typename T>
   inline vec4_t<T> min(vec4_t<T> v, T f) { return vec4_t<T>(min(v.v[0],f),min(v.v[1],f),min(v.v[2],f),min(v.v[3],f)); }
template <typename T>
   inline vec4_t<T> clamp(vec4_t<T> v, T minVal, T maxVal) { return vec4_t<T>(clamp(v.v[0],minVal,maxVal),clamp(v.v[1],minVal,maxVal),clamp(v.v[2],minVal,maxVal),clamp(v.v[3],minVal,maxVal)); }
template <typename T>
   inline vec4_t<T> clamp(vec4_t<T> v, vec4_t<T> minVal, vec4_t<T> maxVal) {
	return vec4_t<T>(clamp(v.v[0],minVal.v[0],maxVal.v[0]),clamp(v.v[1],minVal.v[1],maxVal.v[1]),clamp(v.v[2],minVal.v[2],maxVal.v[2]),clamp(v.v[3],minVal.v[3],maxVal.v[3]));
}
template <typename T>
   inline vec4_t<T> mix(vec4_t<T> x, vec4_t<T> y, float a) { return vec4_t<T>(mix(x.v[0],y.v[0],a),mix(x.v[1],y.v[1],a),mix(x.v[2],y.v[2],a),mix(x.v[3],y.v[3],a)); }
template <typename T>
   inline vec4_t<T> mix(vec4_t<T> x, vec4_t<T> y, vec4_t<T> a) { return vec4_t<T>(mix(x.v[0],y.v[0],a.v[0]),mix(x.v[1],y.v[1],a.v[1]),mix(x.v[2],y.v[2],a.v[2]),mix(x.v[3],y.v[3],a.v[3])); }
template <typename T>
   inline vec4_t<T> step(vec4_t<T> edge, vec4_t<T> x) { return vec4_t<T>(step(edge.v[0],x.v[0]),step(edge.v[1],x.v[1]),step(edge.v[2],x.v[2]),step(edge.v[3],x.v[3])); }
template <typename T>
   inline vec4_t<T> smoothstep(vec4_t<T> edge0,vec4_t<T> edge1,vec4_t<T> x) {
	return vec4_t<T>(smoothstep(edge0.v[0],edge1.v[0],x.v[0]),smoothstep(edge0.v[1],edge1.v[1],x.v[1]),smoothstep(edge0.v[2],edge1.v[2],x.v[2]),smoothstep(edge0.v[3],edge1.v[3],x.v[3]));
}

template <typename T>
   inline T length(vec4_t<T> v) { return static_cast<T>(sqrtf(static_cast<float>(dot(v,v)))); }
template <typename T>
   inline T distance(vec4_t<T> p0, vec4_t<T> p1) { return length(p0-p1); }
template <typename T>
   inline T dot(vec4_t<T> a, vec4_t<T> b) { return a.v[0]*b.v[0] + a.v[1]*b.v[1] + a.v[2]*b.v[2] + a.v[3]*b.v[3]; }
template <typename T>
   inline vec4_t<T> normalize(vec4_t<T> v) { T r = length(v); return vec4_t<T>(v.v[0]/r, v.v[1]/r, v.v[2]/r, v.v[3]/r); }
template <typename T>
   inline vec4_t<T> faceforward(vec4_t<T> N, vec4_t<T> I, vec4_t<T> Nref) {	return dot(Nref,I) < 0 ? N : -N; }
template <typename T>
   inline vec4_t<T> reflect(vec4_t<T> I, vec4_t<T> N) { return I - 2*dot(N,I)*N; }

template <typename T>
   inline bvec4 lessThan(vec4_t<T> a, vec4_t<T> b) { return bvec4(a.v[0]<b.v[0],a.v[1]<b.v[1],a.v[2]<b.v[2],a.v[3]<b.v[3]); }
template <typename T>
   inline bvec4 lessThanEqual(vec4_t<T> a, vec4_t<T> b) { return bvec4(a.v[0]<=b.v[0],a.v[1]<=b.v[1],a.v[2]<=b.v[2],a.v[3]<=b.v[3]); }
template <typename T>
   inline bvec4 greaterThan(vec4_t<T> a, vec4_t<T> b) { return bvec4(a.v[0]>b.v[0],a.v[1]>b.v[1],a.v[2]>b.v[2],a.v[3]>b.v[3]); }
template <typename T>
   inline bvec4 greaterThanEqual(vec4_t<T> a, vec4_t<T> b) { return bvec4(a.v[0]>=b.v[0],a.v[1]>=b.v[1],a.v[2]>=b.v[2],a.v[3]>=b.v[3]); }
template <typename T>
   inline bvec4 equal(vec4_t<T> a, vec4_t<T> b) { return bvec4(a.v[0]==b.v[2],a.v[1]==b.v[1],a.v[2]==b.v[2],a.v[3]==b.v[3]); }
template <typename T>
   inline bvec4 notEqual(vec4_t<T> a, vec4_t<T> b) { return bvec4(a.v[0]!=b.v[2],a.v[1]!=b.v[1],a.v[2]!=b.v[2],a.v[3]!=b.v[3]); }

inline bool any(bvec4 v) { return v.v[0] || v.v[1] || v.v[2] || v.v[3]; }
inline bool all(bvec4 v) { return v.v[0] && v.v[1] && v.v[2] && v.v[3]; }
#ifndef __GNUC__
inline bvec4 not(bvec4 v) { return bvec4(!v.v[0],!v.v[1],!v.v[2],!v.v[3]); }
#endif


/*
 * mat2_t<T> -- 2x2 matrix
 */


template<typename T>
   struct mat2_t
{
	T m[2*2]; // 0..3

	inline mat2_t() { m[0] = m[1] = m[2] = m[3] = 0; }
	inline mat2_t(T f) { m[0] = m[3] = f; m[1] = m[2] = 0; }
	inline mat2_t(T m0, T m1, T m2, T m3) { m[0]=m0; m[1]=m1; m[2]=m2; m[3] = m3; }
};

typedef struct mat2_t<float> mat2;

// operators
template<typename T>
   inline mat2_t<T> operator+(mat2_t<T> a, mat2_t<T> b) { return mat2_t<T>(a.m[0]+b.m[0],a.m[1]+b.m[1],a.m[2]+b.m[2],a.m[3]+b.m[3]); }
template<typename T>
   inline mat2_t<T> operator-(mat2_t<T> a, mat2_t<T> b) { return mat2_t<T>(a.m[0]-b.m[0],a.m[1]-b.m[1],a.m[2]-b.m[2],a.m[3]-b.m[3]); }

template<typename T>
   inline mat2_t<T> operator+(T a, mat2_t<T> m) { return mat2_t<T>(a+m.m[0],a+m.m[1],a+m.m[2],a+m.m[3]); }
template<typename T>
   inline mat2_t<T> operator+(mat2_t<T> m, T a) { return a+m; }

template<typename T>
   inline mat2_t<T> operator-(T a, mat2_t<T> m) { return mat2_t<T>(a-m.m[0],a-m.m[1],a-m.m[2],a-m.m[3]); }
template<typename T>	
   inline mat2_t<T> operator-(mat2_t<T> m, T a) { return mat2_t<T>(m.m[0]-1,m.m[1]-a,m.m[2]-a,m.m[3]-a); }

template<typename T>
   inline mat2_t<T> operator*(T a, mat2_t<T> m) { return mat2_t<T>(a*m.m[0],a*m.m[1],a*m.m[2],a*m.m[3]); }
template<typename T>
   inline mat2_t<T> operator*(mat2_t<T> m, T a) { return a*m; }

template<typename T>
inline vec2_t<T> operator*(mat2_t<T> m, vec2_t<T> v) {
	return vec2_t<T>( 
		m.m[0]*v.v[0] + m.m[2]*v.v[1],
		m.m[1]*v.v[0] + m.m[3]*v.v[1]
	);
}
template<typename T>
inline mat2_t<T> operator*(mat2_t<T> a, mat2_t<T> b) {
	return mat2_t<T>(
		a.m[0]*b.m[0] + a.m[2]*b.m[1],
		a.m[1]*b.m[0] + a.m[3]*b.m[1],
		a.m[0]*b.m[2] + a.m[2]*b.m[3],
		a.m[1]*b.m[2] + a.m[3]*b.m[3]
	);
}

// built-in functions
template<typename T>
inline mat2_t<T> matrixcompmult(mat2_t<T> a, mat2_t<T> b) {
	return mat2_t<T>(a.m[0]*b.m[0],a.m[1]*b.m[1],a.m[2]*b.m[2],a.m[3]*b.m[3]);
}

// additional non-standard functions
template <typename T>
   inline mat2_t<T> transpose(mat2_t<T> m) { return mat2_t<T>(m.m[0],m.m[2],m.m[1],m.m[3]); }

template <typename T>
   inline mat2_t<T> diag(T x, T y) { return mat2_t<T>(x,0,0,y); }
template <typename T>
   inline mat2_t<T> scale(T x, T y) { return mat2_t<T>(x,0,0,y); }


/*
 * mat3_t<T> -- 3x3 matrix
 */


template<typename T>
   struct mat3_t
{
	T m[3*3]; // 0..8

	inline mat3_t() { m[0] = m[1] = m[2] = m[3] = m[4] = m[5] = m[6] = m[7] = m[8] = 0; }
	inline mat3_t(T f) { m[0] = m[4] = m[8] = f; m[1] = m[2] = m[3] = m[5] = m[6] = m[7] = 0; }
	inline mat3_t(T m0,T m1,T m2,T m3,T m4,T m5,T m6,T m7,T m8)	{ m[0]=m0; m[1]=m1; m[2]=m2; m[3]=m3; m[4]=m4; m[5]=m5; m[6]=m6; m[7]=m7; m[8]=m8; }
};

typedef struct mat3_t<float> mat3;

// operators
template<typename T>
   inline mat3_t<T> operator+(mat3_t<T> a, mat3_t<T> b) {
	return mat3_t<T>(a.m[0]+b.m[0], a.m[1]+b.m[1], a.m[2]+b.m[2], a.m[3]+b.m[3], a.m[4]+b.m[4], a.m[5]+b.m[5], a.m[6]+b.m[6], a.m[7]+b.m[7], a.m[8]+b.m[8]);
}
template<typename T>
   inline mat3_t<T> operator-(mat3_t<T> a, mat3_t<T> b) {
	return mat3_t<T>(a.m[0]-b.m[0], a.m[1]-b.m[1], a.m[2]-b.m[2], a.m[3]-b.m[3], a.m[4]-b.m[4], a.m[5]-b.m[5], a.m[6]-b.m[6], a.m[7]-b.m[7], a.m[8]-b.m[8]);
}

template<typename T>
   inline mat3_t<T> operator+(T a, mat3_t<T> m) { return mat3_t<T>(a+m.m[0], a+m.m[1], a+m.m[2], a+m.m[3], a+m.m[4], a+m.m[5], a+m.m[6], a+m.m[7], a+m.m[8]); }
template<typename T>
   inline mat3_t<T> operator+(mat3_t<T> m, T a) { return a+m; }

template <typename T>
   inline mat3_t<T> operator-(T a, mat3_t<T> m) { return mat3_t<T>(a-m.m[0], a-m.m[1], a-m.m[2], a-m.m[3], a-m.m[4], a-m.m[5], a-m.m[6], a-m.m[7], a-m.m[8]); }
template <typename T>
   inline mat3_t<T> operator-(mat3_t<T> m, T a) { return mat3_t<T>(m.m[0]-a, m.m[1]-a, m.m[2]-a, m.m[3]-a, m.m[4]-a, m.m[5]-a, m.m[6]-a, m.m[7]-a, m.m[8]-a); }

template <typename T>
   inline mat3_t<T> operator*(T a, mat3_t<T> m) { return mat3_t<T>(a*m.m[0], a*m.m[1], a*m.m[2], a*m.m[3], a*m.m[4], a*m.m[5], a*m.m[6], a*m.m[7], a*m.m[8]); }
template <typename T>
   inline mat3_t<T> operator*(mat3_t<T> m, T a) { return a*m; }

template <typename T>
   inline vec3_t<T> operator*(mat3_t<T> m, vec3_t<T> v) {
	return vec3_t<T>(
					 m.m[0]*v.v[0] + m.m[3]*v.v[1] + m.m[6]*v.v[2],
					 m.m[1]*v.v[0] + m.m[4]*v.v[1] + m.m[7]*v.v[2],
					 m.m[2]*v.v[0] + m.m[5]*v.v[1] + m.m[8]*v.v[2]
					);
}
template <typename T>
   inline mat3_t<T> operator*(mat3_t<T> a, mat3_t<T> b) {
	return mat3_t<T>(
					 a.m[0]*b.m[0] + a.m[3]*b.m[1] + a.m[6]*b.m[2],
					 a.m[1]*b.m[0] + a.m[4]*b.m[1] + a.m[7]*b.m[2],
					 a.m[2]*b.m[0] + a.m[5]*b.m[1] + a.m[8]*b.m[2],

					 a.m[0]*b.m[3] + a.m[3]*b.m[4] + a.m[6]*b.m[5],
					 a.m[1]*b.m[3] + a.m[4]*b.m[4] + a.m[7]*b.m[5],
					 a.m[2]*b.m[3] + a.m[5]*b.m[4] + a.m[8]*b.m[5],

					 a.m[0]*b.m[6] + a.m[3]*b.m[7] + a.m[6]*b.m[8],
					 a.m[1]*b.m[6] + a.m[4]*b.m[7] + a.m[7]*b.m[8],
					 a.m[2]*b.m[6] + a.m[5]*b.m[7] + a.m[8]*b.m[8] 
					);
}

// built-in functions
template <typename T>
   inline mat3_t<T> matrixcompmult(mat3_t<T> a, mat3_t<T> b) {
	return mat3_t<T>(a.m[0]*b.m[0], a.m[1]*b.m[1], a.m[2]*b.m[2], a.m[3]*b.m[3], a.m[4]*b.m[4], a.m[5]*b.m[5], a.m[6]*b.m[6], a.m[7]*b.m[7], a.m[8]*b.m[8]);
}

// additional non-standard functions
template <typename T>
   inline mat3_t<T> transpose(mat3_t<T> m) {
	return mat3_t<T>(m.m[0], m.m[3], m.m[6], m.m[1], m.m[4], m.m[7], m.m[2], m.m[5], m.m[8]);
}
template <typename T>
   inline mat3_t<T> diag(T x,T y,T z) { return mat3_t<T>(x,0,0,0,y,0,0,0,z); }
template <typename T>
   inline mat3_t<T> scale(T x,T y,T z) { return mat3_t<T>(x,0,0,0,y,0,0,0,z); }



/*
 * mat4_t<T> --- 4x4 matrix
 */



template<typename T>
   struct mat4_t
{
	T m[4*4]; // 0..15

	inline mat4_t()    { m[0] = m[1] = m[2] = m[3] = m[4] = m[5] = m[6] = m[7] = m[8] = m[9] = m[10] = m[11] = m[12] = m[13] = m[14] = m[15] = 0; }
	inline mat4_t(T f) { m[0] = m[5] = m[10] = m[15] = f; m[1] = m[2] = m[3] = m[4] = m[6] = m[7] = m[8] = m[9] = m[11] = m[12] = m[13] = m[14] = 0; }
	inline mat4_t(T m0, T m1, T m2, T m3, T m4, T m5, T m6, T m7, T m8, T m9, T m10, T m11, T m12, T m13, T m14, T m15) {
		m[0]=m0; m[1]=m1; m[2]=m2; m[3]=m3; m[4]=m4; m[5]=m5; m[6]=m6; m[7]=m7; m[8]=m8; m[9]=m9; m[10]=m10; m[11]=m11; m[12]=m12; m[13]=m13; m[14]=m14; m[15]=m15;
	}
};

typedef struct mat4_t<float> mat4;

// operators
template <typename T>
   inline mat4_t<T> operator+(mat4_t<T> a, mat4_t<T> b) {
	return mat4_t<T>(a.m[0]+b.m[0], a.m[1]+b.m[1], a.m[2]+b.m[2], a.m[3]+b.m[3], a.m[4]+b.m[4], a.m[5]+b.m[5], a.m[6]+b.m[6], a.m[7]+b.m[7], a.m[8]+b.m[8], a.m[9]+b.m[9],
					 a.m[10]+b.m[10], a.m[11]+b.m[11], a.m[12]+b.m[12], a.m[13]+b.m[13], a.m[14]+b.m[14], a.m[15]+b.m[15]);
}
template <typename T>
   inline mat4_t<T> operator-(mat4_t<T> a, mat4_t<T> b) {
	return mat4_t<T>(a.m[0]-b.m[0], a.m[1]-b.m[1], a.m[2]-b.m[2], a.m[3]-b.m[3], a.m[4]-b.m[4], a.m[5]-b.m[5], a.m[6]-b.m[6], a.m[7]-b.m[7], a.m[8]-b.m[8], a.m[9]-b.m[9],
					 a.m[10]-b.m[10], a.m[11]-b.m[11], a.m[12]-b.m[12], a.m[13]-b.m[13], a.m[14]-b.m[14], a.m[15]-b.m[15]);
}

template <typename T>
   inline mat4_t<T> operator+(T a, mat4_t<T> m) {
	return mat4_t<T>(a+m.m[0],a+m.m[1],a+m.m[2],a+m.m[3],a+m.m[4],a+m.m[5],a+m.m[6],a+m.m[7],a+m.m[8],a+m.m[9],a+m.m[10],a+m.m[11],a+m.m[12],a+m.m[13],a+m.m[14],a+m.m[15]);
}
template <typename T>
   inline mat4_t<T> operator+(mat4_t<T> m, T a) { return a+m; }

template <typename T>
   inline mat4_t<T> operator-(T a, mat4_t<T> m) {
	return mat4_t<T>(a-m.m[0],a-m.m[1],a-m.m[2],a-m.m[3],a-m.m[4],a-m.m[5],a-m.m[6],a-m.m[7],a-m.m[8],a-m.m[9],a-m.m[10],a-m.m[11],a-m.m[12],a-m.m[13],a-m.m[14],a-m.m[15]);
}
template <typename T>
   inline mat4_t<T> operator-(mat4_t<T> m, T a) {
	return mat4_t<T>(m.m[0]-a,m.m[1]-a,m.m[2]-a,m.m[3]-a,m.m[4]-a,m.m[5]-a,m.m[6]-a,m.m[7]-a,m.m[8]-a,m.m[9]-a,m.m[10]-a,m.m[11]-a,m.m[12]-a,m.m[13]-a,m.m[14]-a,m.m[15]-a);
}

template <typename T>
   inline mat4_t<T> operator*(T a, mat4_t<T> m) {
	return mat4_t<T>(a*m.m[0],a*m.m[1],a*m.m[2],a*m.m[3],a*m.m[4],a*m.m[5],a*m.m[6],a*m.m[7],a*m.m[8],a*m.m[9],a*m.m[10],a*m.m[11],a*m.m[12],a*m.m[13],a*m.m[14],a*m.m[15]);
}
template <typename T>
   inline mat4_t<T> operator*(mat4_t<T> m, T a) { return a*m; }

template <typename T>
   inline vec4_t<T> operator*(mat4_t<T> m, vec4_t<T> v) {
	return vec4_t<T>(
					 m.m[0]*v.v[0] + m.m[4]*v.v[1] + m.m[8]*v.v[2] + m.m[12]*v.v[3],
					 m.m[1]*v.v[0] + m.m[5]*v.v[1] + m.m[9]*v.v[2] + m.m[13]*v.v[3],
					 m.m[2]*v.v[0] + m.m[6]*v.v[1] + m.m[10]*v.v[2]+ m.m[14]*v.v[3],
					 m.m[3]*v.v[0] + m.m[7]*v.v[1] + m.m[11]*v.v[2]+ m.m[15]*v.v[3]
					);
}
template <typename T>
   inline mat4_t<T> operator*(mat4_t<T> a, mat4_t<T> b) {
	return mat4_t<T>(
					 a.m[0]*b.m[0] + a.m[4]*b.m[1] + a.m[8]*b.m[2] + a.m[12]*b.m[3],
					 a.m[1]*b.m[0] + a.m[5]*b.m[1] + a.m[9]*b.m[2] + a.m[13]*b.m[3],
					 a.m[2]*b.m[0] + a.m[6]*b.m[1] + a.m[10]*b.m[2]+ a.m[14]*b.m[3],
					 a.m[3]*b.m[0] + a.m[7]*b.m[1] + a.m[11]*b.m[2]+ a.m[15]*b.m[3],

					 a.m[0]*b.m[4] + a.m[4]*b.m[5] + a.m[8]*b.m[6] + a.m[12]*b.m[7],
					 a.m[1]*b.m[4] + a.m[5]*b.m[5] + a.m[9]*b.m[6] + a.m[13]*b.m[7],
					 a.m[2]*b.m[4] + a.m[6]*b.m[5] + a.m[10]*b.m[6]+ a.m[14]*b.m[7],
					 a.m[3]*b.m[4] + a.m[7]*b.m[5] + a.m[11]*b.m[6]+ a.m[15]*b.m[7],

					 a.m[0]*b.m[8] + a.m[4]*b.m[9] + a.m[8]*b.m[10] + a.m[12]*b.m[11],
					 a.m[1]*b.m[8] + a.m[5]*b.m[9] + a.m[9]*b.m[10] + a.m[13]*b.m[11],
					 a.m[2]*b.m[8] + a.m[6]*b.m[9] + a.m[10]*b.m[10]+ a.m[14]*b.m[11],
					 a.m[3]*b.m[8] + a.m[7]*b.m[9] + a.m[11]*b.m[10]+ a.m[15]*b.m[11],

					 a.m[0]*b.m[12] + a.m[4]*b.m[13] + a.m[8]*b.m[14] + a.m[12]*b.m[15],
					 a.m[1]*b.m[12] + a.m[5]*b.m[13] + a.m[9]*b.m[14] + a.m[13]*b.m[15],
					 a.m[2]*b.m[12] + a.m[6]*b.m[13] + a.m[10]*b.m[14]+ a.m[14]*b.m[15],
					 a.m[3]*b.m[12] + a.m[7]*b.m[13] + a.m[11]*b.m[14]+ a.m[15]*b.m[15] 

					);
}

// built-in functions
template <typename T>
   inline mat4_t<T> matrixcompmult(mat4_t<T> a, mat4_t<T> b) {
	return mat4_t<T>(a.m[0]*b.m[0],a.m[1]*b.m[1],a.m[2]*b.m[2],a.m[3]*b.m[3],a.m[4]*b.m[4],a.m[5]*b.m[5],a.m[6]*b.m[6],a.m[7]*b.m[7],a.m[8]*b.m[8],
					 a.m[9]*b.m[9],a.m[10]*b.m[10],a.m[11]*b.m[11],a.m[12]*b.m[12],a.m[13]*b.m[13],a.m[14]*b.m[14],a.m[15]*b.m[15]);
}

// additional non-standard functions
template <typename T>
   inline mat4_t<T> transpose(mat4_t<T> m) {
	return mat4_t<T>(m.m[0],m.m[4],m.m[8],m.m[12],m.m[1],m.m[5],m.m[9],m.m[13],m.m[2],m.m[6],m.m[10],m.m[14],m.m[3],m.m[7],m.m[11],m.m[15]);
}
template <typename T>
   inline mat4_t<T> diag(T x,T y, T z, T w) {
	return mat4_t<T>(x,0,0,0, 0,y,0,0, 0,0,z,0, 0,0,0,w);
}
template <typename T>
   inline mat4_t<T> scale(T x,T y, T z, T w) {
	return mat4_t<T>(x,0,0,0, 0,y,0,0, 0,0,z,0, 0,0,0,w);
}


/*
 * dimension downgrade conversions
 */


template <typename T> inline vec2_t<T>::vec2_t(vec3_t<T> v) { v[0]=v.v[0]; v[1]=v.v[1]; }
template <typename T> inline vec3_t<T>::vec3_t(vec4_t<T> v) { v[0]=v.v[0]; v[1]=v.v[1]; v[2]=v.v[2]; }


// init_xxx() funcs for compatibility with glsl_gen.h
template <typename T> inline vec2_t<T> init_vec2(T v0, T v1) { return vec2_t<T>(v0,v1); }
template <typename T> inline vec3_t<T> init_vec3(T v0, T v1, T v2) { return vec3_t<T>(v0,v1,v2); }
template <typename T> inline vec4_t<T> init_vec4(T v0, T v1, T v2, T v3) { return vec4_t<T>(v0,v1,v2,v3); }

template <typename T> inline mat2_t<T> init_mat2(T m0, T m1, T m2, T m3) { 
	return mat2_t<T>(m0,m1,m2,m3);
}
template <typename T> inline mat3_t<T> init_mat3(T m0,T m1,T m2,T m3,T m4,T m5,T m6,T m7,T m8) { 
	return mat3_t<T>(m0,m1,m2,m3,m4,m5,m6,m7,m8);
}
template <typename T> inline mat3_t<T> init_mat4(T m0, T m1, T m2, T m3, T m4, T m5, T m6, T m7, T m8, T m9, T m10, T m11, T m12, T m13, T m14, T m15) { 
	return mat3_t<T>(m0,m1,m2,m3,m4,m5,m6,m7,m8,m9,m10,m11,m12,m13,m14,m15);
}

#endif /* INCLUDED_GLSL_H */

